Daily Exercise 25

Ecosystem Science and Sustainability 330

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.5.1     ✔ tibble    3.2.1
✔ lubridate 1.9.4     ✔ tidyr     1.3.1
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(sf)

Warning: package 'sf' was built under R version 4.4.3

Linking to GEOS 3.13.0, GDAL 3.10.1, PROJ 9.5.1; sf_use_s2() is TRUE

library(AOI)
library(tmap)

Warning: package 'tmap' was built under R version 4.4.3

library(tigris)

Warning: package 'tigris' was built under R version 4.4.3

To enable caching of data, set `options(tigris_use_cache = TRUE)`
in your R script or .Rprofile.

Attaching package: 'tigris'

The following object is masked from 'package:AOI':

    list_states

rivers <- read_sf("data/MajorRivers.shp")

rivers

Simple feature collection with 98 features and 4 fields
Geometry type: MULTILINESTRING
Dimension:     XY
Bounding box:  xmin: -164.8874 ymin: -36.96945 xmax: 160.7636 ymax: 71.39249
Geodetic CRS:  WGS 84
# A tibble: 98 × 5
   NAME          SYSTEM MILES KILOMETERS                                geometry
   <chr>         <chr>  <dbl>      <dbl>                   <MULTILINESTRING [°]>
 1 Kolyma        <NA>   2552.      4106. ((144.8419 61.75915, 144.8258 61.8036,…
 2 Parana        Parana 1616.      2601. ((-51.0064 -20.07941, -51.02972 -20.22…
 3 San Francisco <NA>   1494.      2404. ((-46.43639 -20.25807, -46.49835 -20.2…
 4 Japura        Amazon 1223.      1968. ((-76.71056 1.624166, -76.70029 1.6883…
 5 Putumayo      Amazon  890.      1432. ((-76.86806 1.300553, -76.86695 1.295,…
 6 Rio Maranon   Amazon  889.      1431. ((-73.5079 -4.459834, -73.79197 -4.621…
 7 Ucayali       Amazon 1298.      2089. ((-73.5079 -4.459834, -73.51585 -4.506…
 8 Guapore       Amazon  394.       634. ((-65.39585 -10.39333, -65.39578 -10.3…
 9 Madre de Dios Amazon  568.       914. ((-65.39585 -10.39333, -65.45279 -10.4…
10 Amazon        Amazon 1890.      3042. ((-73.5079 -4.459834, -73.45141 -4.427…
# ℹ 88 more rows

mississippi <- rivers %>%
  filter(str_detect(NAME, "Mississippi|Missouri|Ohio|Arkansas|Red River"))

options(tigris_use_cache = TRUE)
counties <- counties(cb = TRUE, year = 2020) %>%
  st_transform(crs = st_crs(rivers))

counties_mississippi <- st_filter(counties, mississippi)

tmap_mode("view")

ℹ tmap mode set to "view".

tm_shape(counties_mississippi) + 
  tm_borders(col = "black") +
tm_shape(mississippi) +
  tm_lines(col = "blue")

library(readr)

cities <- read_csv("data/uscities.csv") 

Rows: 31254 Columns: 17
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (9): city, city_ascii, state_id, state_name, county_fips, county_name, s...
dbl (6): lat, lng, population, density, ranking, id
lgl (2): military, incorporated

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

library(sf)

cities_sf <- st_as_sf(cities, coords = c("lng", "lat"), crs = 4326)  
cities_sf <- st_transform(cities_sf, crs = st_crs(counties))  

intersecting_counties <- counties %>%
  st_filter(mississippi, .predicate = st_intersects)

city_matches <- st_join(cities_sf, intersecting_counties, left = FALSE)

library(dplyr)

county_pop <- city_matches %>%
  group_by(GEOID) %>%  
  summarize(total_pop = sum(population, na.rm = TRUE))

intersecting_counties <- intersecting_counties %>%
  left_join(st_drop_geometry(county_pop), by = "GEOID")

library(ggplot2)

ggplot() +
  geom_sf(data = mississippi, color = "blue", size = 0.5) +
  
  geom_sf(data = intersecting_counties, aes(fill = total_pop), color = "black", size = 0.1) +
  
  scale_fill_viridis_c(option = "C", name = "Total Urban Population") +
  theme_minimal() +
  labs(title = "Counties Intersecting the Mississippi River System",
       subtitle = "Color-coded by Total Urban Population") +
  theme(legend.position = "bottom")